Optosurf simulation

Loading BokehJS ...
Code
def subplot3d(files, offaxis, onaxis, rows, cols, path):
    """
    Generates a 3D plot of the data

    Parameters
    ----------
    files (list): List of files to plot
    df_slice (dataframe): Dataframe with the slices to plot
    
    Returns
    -------
    subplot (plotly.graph_objects.Figure): Plotly figure with the 3D plot
    spots (dict): Dictionary with the data of the files
    """
    # a. Define on and off axis as well as subplots
    subplot = make_subplots(rows=rows, cols=cols, subplot_titles=files, 
                    specs=[[{'type': 'surface'}]*cols]*rows)
    x = np.linspace(0, 10, 100)
    z = np.linspace(0, 10, 100)
    y = np.linspace(-15.5, 15.5, 100)
    X, Z = np.meshgrid(x, z)
    Y, Z2 = np.meshgrid(y, z)
    count = 0
    spots = {}
    col = 1

    # b. Read file and create matrix
    for i, file in enumerate(files):
        file_path = path + file
        # c. Read the .csv file
        vals = np.genfromtxt(file_path, delimiter=',')
        
        # d. Normalize values
        normalized_vals = vals[:, 1:] 
        # normalized_vals = 10 * vals[:, 1:] / np.max(np.max(vals[:, 1:]))
        
        # e. 3D surface plot
        surface = go.Surface(x=offaxis, y=onaxis, z=normalized_vals, colorscale='jet', 
                   showscale=False)
        row = i // 3 + 1  # Compute the row index based on the file index
        col = i % 3 + 1   # Compute the column index based on the file index
        subplot.add_trace(surface, row=row, col=col)
        
        # f. Create dictionary
        spots[file] = vals
    return subplot, spots

file1 = 'refrockley_focus_0.csv'
file2 = 'refrockley_focus_m1.csv' 
file3 = 'refrockley_focus_m2.csv'
file4 = 'refrockley_focus_m3.csv'
file5 = 'refrockley_focus_p1.csv'
file6 = 'refrockley_focus_p2.csv'
files = [file1, file2, file3, file4, file5, file6]
offaxis = np.arange(-2048, 2048, 1)*24/4096 + 1.06
onaxis=np.arange(-15.5,16.5,1)

path = '/Users/carlosreyes/Desktop/Academics/postdoc/confirm/python/optosurf/Streamlit/data/f/refocus/'
subplot, spots = subplot3d([file1, file2, file3, file4, file5, file6], offaxis, onaxis, 2, 3, path)
subplot.update_scenes(camera_eye=dict(x=-0.75, y=-1.0, z=1.25))
subplot.update_layout(width=1200, height=900)
subplot
slice1 = 100
slice2 = 4000
print(offaxis[slice1], offaxis[slice2])
-10.3540625 12.4975

1 Datasets

2 1st minimization iteration with original optosurf axis

Code
# 5. Define cost function
def cost_function(params, y, x):
    """
    Generates a 3D plot of the data

    Parameters
    ----------
    params (list): List of parameters to fit
    y (array): Array with the data to fit
    
    Returns
    -------
    rmse (float): Root mean square error of the fit0.
    """
    x0, A0 = params
    x_new = x + x0
    y_modified = A0*pchip(x_new)

    mse = np.mean((y - y_modified) ** 2)
    rmse = np.sqrt(mse)
    return rmse

# slices = np.arange(20, 4090, step=1)
slices = np.arange(slice1, slice2+1, step=1)

# methods = ['Powell', 'CG', 'L-BFGS-B', 'SLSQP', 'trust-constr']
method = 'Powell'

# Get initial base function
smooth_df = pd.read_csv("data/f/smooth_df.csv")
x_base = smooth_df['xaxis']
y_base = smooth_df['yaxis']
pchip = PchipInterpolator(x_base, y_base)

offaxis = np.arange(-2048, 2048, 1)*24/4096 + 1.06
onaxis=np.arange(-15.5,16.5,1)

minimized_df_dict = {}
files = [file1, file2, file3, file4, file5, file6]

for file in files:
    minimized_df = pd.DataFrame(columns=["angle", "x0", "difference", "slice", "amplitude", "rmse"])
    dataset = spots[file]
    for j, slice in enumerate(slices):
        x0 = offaxis[slice]
        Abase = 1.0
        guess_no_back = [x0, Abase]

        # Call minimization function
        y = dataset[:,slice]
        cost_fn = lambda p:cost_function(p, y, onaxis)

        result = minimize(cost_fn, guess_no_back, method=method,)
        optimized_parameters = list(result.x)
        x0_opt, A0_opt, = optimized_parameters

        # Calculate optimized function
        x_new_opt = onaxis + (x0_opt)
        y_optimized = A0_opt*pchip(x_new_opt)

        # Calculate error
        mse = np.mean((y - y_optimized) ** 2)
        rmse = np.sqrt(mse)
        row = [offaxis[slice], x0_opt, offaxis[slice] - x0_opt, slice, A0_opt, rmse]
        minimized_df.loc[j] = row
    minimized_df_dict[file] = minimized_df.copy()
Code
sl = 20
yexp = dataset[:,sl]
print(offaxis[sl])
pchipplot = figure(title='pchip plot', width= 700, height= 400, x_range = (-15,15))
pchipplot.line(onaxis, yexp, line_width=2, color=new_colors[0], legend_label = f'Experimental data ({sl})')
pchipplot.circle(onaxis, yexp, size  = 5, legend_label = f'Experimental data ({sl})', color = new_colors[0])

pchipplot.line(x_base, y_base, line_width=2, color=new_colors[1], legend_label = f'Base function')
pchipplot.circle(x_base, y_base, size  = 5, legend_label = f'Base function', color= new_colors[1])

yonaxis = pchip(onaxis)
pchipplot.line(onaxis, yonaxis, line_width=2, color=new_colors[2], legend_label = 'pchip base')
pchipplot.circle(onaxis, yonaxis, size  = 5, legend_label = 'pchip base', color=new_colors[2])

x_new = onaxis - 10
y_modified = 0.95*pchip(x_new)
pchipplot.line(onaxis, y_modified, line_width=2, color=new_colors[3], dash = 'dashed', legend_label = 'x0 = 15')
pchipplot.circle(onaxis, y_modified, size  = 5, color = new_colors[3], legend_label = 'x0 = 15')

pchipplot = plot_format(pchipplot, "angle", "amplitude", "top_left", "10pt", "11pt", "8pt")
show(pchipplot)

mse = np.mean((yexp - y_modified) ** 2)
rmse = np.sqrt(mse)
print(rmse)
-10.8228125
2280.254029991687

2.1 Shifted axis and amplitude correction factor

Code
p1 = figure(title=f'a. Rotation angle vs. x0')
p2 = figure(title='c. Rotation angle vs. RMSE')
p3 = figure(title='b. x0 vs. difference (mechanical angle - x0)')
p4 = figure(title='d. Rotation angle vs. amplitude')
p5 = figure(title='Shifted optosurf axis', width=1000, height=350, x_range=(-28,28))

def linear_function(x, m , b):
    return m*x + b
shifted_axis_dict = {}
pchip_amps = {}
angles_interp_dict = {}
amplitudes_interp_dict = {}

for z, file in enumerate(files):
    # a. Rotation angle vs. x0
    p1.line(x=minimized_df_dict[file]['angle'], y=minimized_df_dict[file]['x0'], line_width=2, color = new_colors[z], legend_label=file)

    # b. Rotation angle vs RMSE
    p2.line(x=minimized_df_dict[file]['angle'], y=minimized_df_dict[file]['rmse'], line_width=2, color = new_colors[z], legend_label = file)

    # c. Polynomial fit of x0 vs difference
    x0 = minimized_df_dict[file]['x0'].values
    angles = minimized_df_dict[file]['angle'].values
    difference = minimized_df_dict[file]['difference'].values
    poly = PolynomialFeatures(degree=3)
    # x0_poly = poly.fit_transform(angles.reshape(-1,1))
    x0_poly = poly.fit_transform(x0.reshape(-1,1))
    poly_model = LinearRegression().fit(X=x0_poly, y=difference)
    xaxis = np.arange(-15.5, 15.5, 0.1)
    xaxis_poly = poly.transform(xaxis.reshape(-1,1))
    ypoly = poly_model.predict(xaxis_poly)
    p3.line(x=x0, y=difference, line_width=2, color = new_colors[z], legend_label = file)
    # p3.line(x=xaxis, y=ypoly, line_width=2, color = new_colors[4], dash='dashed')
    
    # d. Calculate shifted axis
    onaxis_poly = poly.transform(onaxis.reshape(-1,1))
    y_shifted = poly_model.predict(onaxis_poly)
    p3.circle(onaxis, y_shifted, size = 6, fill_color = new_colors[z], color='blue')
    shifted_axis = onaxis + y_shifted
    shifted_axis_dict[file] = shifted_axis
    
    # e. Rotation angle vs amplitude with interpolation
    angles = minimized_df_dict[file]['angle'].values
    amplitudes = minimized_df_dict[file]['amplitude'].values
    p4.line(x=angles, y=amplitudes, line_width=1, color = new_colors[z], legend_label = file)
    
    angles = np.insert(angles, 0, -20)
    angles = np.append(angles, 20)
    amplitudes = np.insert(amplitudes, 0, amplitudes[0])
    amplitudes = np.append(amplitudes, amplitudes[-1])
    pchip_amp = PchipInterpolator(angles, amplitudes)
    pchip_amps[file] = pchip_amp
    angles_interp = np.arange(-18, 18, 0.0005)
    # amplitude_interp = pchip_amp(angles_interp)
    amplitude_interp = pchip_amp(shifted_axis)
    angles_interp_dict[file] = angles_interp
    amplitudes_interp_dict[file] = amplitude_interp
    # p4.line(x=angles_interp, y=amplitude_interp, line_width=2, color = new_colors[2], dash='dashed')
    p4.circle(x=shifted_axis, y=amplitude_interp, size=5, color = new_colors[2],)

    p5.circle(x=shifted_axis, y=np.zeros(len(onaxis))+0.5*z, line_width=2, color = new_colors[z], size = 7, legend_label=file)

# f. Plots format
p1.xaxis.ticker.desired_num_ticks = 10
p1 = plot_format(p1, "Rotation angle", "x0", "top_left", "10pt", "8pt", "8pt")
p2.xaxis.ticker.desired_num_ticks = 10
p2 = plot_format(p2, "Rotation angle", "RMSE", "top_left", "10pt", "11pt", "8pt")
p3.xaxis.ticker.desired_num_ticks = 10
p3 = plot_format(p3, "x0", "difference", "top_left", "10pt", "11pt", "8pt")
p4.xaxis.ticker.desired_num_ticks = 10
p4 = plot_format(p4, "Rotation angle", "amplitude", "bottom_left", "10pt", "11pt", "8pt")

p5.circle(x=onaxis, y=np.zeros(len(onaxis))-0.5, line_width=2, color = new_colors[z+1], size = 7, legend_label='Original axis')
p5 = plot_format(p5, "Angle", "", "top_left", "10pt", "11pt", "8pt")
p5.xaxis.ticker.desired_num_ticks = 20

grid = gridplot(children=[p1, p2, p3, p4], ncols=2, merge_tools=False, width = 750, height = 400)
show(grid)
Code
show(p5)

2.2 Base function recreation

Code
# g. Reconstruct base function
# base_1 = 20
# base_2 = 4088
base_1 = slice1
base_2 = slice2
data_str = 'refrockley_focus_0.csv' # 'refrockley_focus_0.csv', 'refrockley_focus_m3.csv' , 'refrockley_focus_p2.csv'
base_step = 1
minimized_df_2 = minimized_df_dict[data_str].copy()
shifted_axis = shifted_axis_dict[data_str]
pchip_amp = pchip_amps[data_str]

slices_base = np.arange(base_1, base_2+1, step=base_step)
minimized_df_2 = minimized_df_2.set_index(['slice'])

p6 = figure(title= 'Modified slices', width = 1000, height = 350,)
x_base = []
y_base = []
slice_array = []
vals = spots[data_str]
experimental_data = spots[data_str]
vals_it_0 = np.zeros_like(vals)

amp = pchip_amp(shifted_axis)
color = 0 
for k, slice in enumerate(slices_base):
    # Get the starting experimental 32 points
    y_temp = vals[:,slice]
    angle = minimized_df_2.loc[slice, 'angle']
    x0 = minimized_df_2.loc[slice, 'x0']
    # yamp= y_temp/amp
    yamp= y_temp*amp
    # yamp = y_temp
    if slice % 250 == 0:
        p6.circle(x=shifted_axis, y=y_temp, size=4, color = new_colors[color], legend_label = f'Original slice {slice}', )
        p6.line(x=shifted_axis, y=y_temp, line_width=2, color = new_colors[color], legend_label = f'Original slice {slice}', line_dash = 'dashed')

        p6.circle(x=shifted_axis, y=yamp, size=4, color = new_colors[color+1], legend_label = f'Modified slice {slice}')
        p6.line(x=shifted_axis, y=yamp, line_width=2, color = new_colors[color+1], legend_label = f'Modified slice {slice}')
        color += 1
        
    # Reshift to the center
    shifted_axis_2 = shifted_axis + angle
    x_base.extend(list(shifted_axis_2))
    y_base.extend(list(yamp))
    slice_array.extend(np.ones(len(shifted_axis_2))*slice)

    # update vals
    # vals[:,slice] = yamp
    vals_it_0[:,slice] = yamp

# vals = vals_new.copy()
# print(np.max(vals_new))
p6 = plot_format(p6, "Angle", "Intensity", "top_right", "10pt", "11pt", "8pt")
show(p6)
Code
# 2. Create base function df
base_function_df = pd.DataFrame({'xaxis': x_base, 'yaxis': y_base, 'slice': slice_array})
base_function_df = base_function_df.sort_values(by='xaxis')

base_function_merged = base_function_df.merge(minimized_df_2, on=['slice', 'slice'])
base_function_merged.drop(columns=['rmse', 'difference'], inplace=True)
base_function_merged['yaxisamp'] = base_function_merged['yaxis'] / base_function_merged['amplitude']

iteration_1['minimized_df'] = minimized_df_dict[data_str].copy()
iteration_1['shifted'] = shifted_axis_dict[data_str].copy()
iteration_1['base'] = base_function_merged.copy()

all_filtered = np.arange(1000, 3000, step=4)
base_section = base_function_merged.copy()
base_section = base_section[base_section['slice'].isin(all_filtered)]
base_section = base_section.sort_values(by=['xaxis'])

p7 = figure(title=f'Base function {offaxis[1000]:.3f} to {offaxis[3000]:.3f}', width=850, height=650, x_range=(-1.5, 1.5), y_range = (25000, 38000))
p7.circle(x=base_section['xaxis'], y=base_section['yaxis'], line_width=2, color = new_colors[0], size = 3, legend_label='Base function')
p7 = plot_format(p7, "Angle", "Amplitude", "top_left", "10pt", "11pt", "8pt")
show(p7)

3 Minimization iteration with original optosurf axis (refrockley_focus_0.csv)

Code
new_colors = []
for i in range(42):
        new_colors.append('#9D6C97')
        new_colors.append('#9DC3E6')
        new_colors.append('#9DD9C5')
        new_colors.append('#F2DDA4')
        new_colors.append('#C4A287')
        new_colors.append('#F49097')
        new_colors.append('#E5625E')
        new_colors.append('#80FFE8')
        new_colors.append('#F7EF81')


def cost_function_a(params, y, x, chip):
    """
    Generates a 3D plot of the data

    Parameters
    ----------
    params (list): List of parameters to fit
    y (array): Array with the data to fit
    
    Returns
    -------
    rmse (float): Root mean square error of the fit0.
    """
    x0, A0 = params
    x_new = x + x0
    y_modified = A0*chip(x_new)
    mse = np.mean((y - y_modified) ** 2)
    rmse = np.sqrt(mse)
    return rmse

def minimize_procedure(slices_list, starting_minimized, x, chip, values):
    """
    Minimizes the cost function for each slice from the 3D dataset

    Parameters
    ----------
    slices (list): List of slices to optimize
    starting_minimized (DataFrame): DataFrame with initial guesses for the parameters

    Returns
    -------
    dfc (DataFrame): DataFrame with the optimized parameters and RMSE for each slice.
    """
    dfc = pd.DataFrame(columns=["angle", "x0", "difference", "slice", "amplitude", "rmse"])
    for j, sliceb in enumerate(slices_list):
        # Get initial guesses from previous df
        x0 = starting_minimized.loc[sliceb]['x0']
        angle = starting_minimized.loc[sliceb]['angle']
        Abase = starting_minimized.loc[sliceb]['amplitude']
        guess = [x0, Abase]

        # Call minimize function
        ya = values[:, sliceb]
        cost_fn = lambda p:cost_function_a(p, ya, x, chip)
        result = minimize(cost_fn, guess, method=method)
        optimized_parameters = list(result.x)
        x0_opt, A0_opt, = optimized_parameters
        
        # Calculate optimized function
        x_new_opt = x + (x0_opt)
        y_optimized = A0_opt*chip(x_new_opt) 

        # Calculate error
        mse = np.mean((ya - y_optimized) ** 2)
        rmse = np.sqrt(mse)
        row = [offaxis[sliceb], x0_opt, offaxis[sliceb] - x0_opt, sliceb, A0_opt, rmse]
        dfc.loc[j] = row
    return dfc

def polyfit(degree, dfa):
    """
    Performs polynomial regression with a given degree on the 'x0' and 'difference' columns of a given dataframe 
    and returns the x-axis, y-predictions, and trained regression model.

    Parameters
    ----------
    degree (int): Degree of the polynomial regression
    dfa (DataFrame): DataFrame with 'x0' and 'difference' columns for the regression

    Returns
    -------
    xaxis (array): Array with x-axis values
    ypredictions (array): Array with y-predictions for each x-axis value
    model (LinearRegression object): Trained linear regression model
    """
    poly = PolynomialFeatures(degree=degree)
    x0 = dfa['angle'].values
    difference = dfa['difference'].values
    x0_poly = poly.fit_transform(x0.reshape(-1,1))
    model = LinearRegression().fit(X=x0_poly, y=difference)
    xaxis = np.arange(-15.5, 15.5, 0.1)
    xaxis_poly = poly.transform(xaxis.reshape(-1,1))
    ypredictions = model.predict(xaxis_poly)
    return xaxis, ypredictions, model


def filter_base(base_function, slice_filter):
    mask = base_function['slice'].isin(slice_filter)
    base_function = base_function.loc[mask].reset_index()
    base_function.loc[0] = [0, -500, 0, 0, 0, 0, 0, 0]
    base_function.loc[-1] = [0, 500, 0, 0, 0, 0, 0, 0]
    base_function.sort_values(by=['xaxis'], inplace=True)
    return base_function


def average_base(base_function, window_size):
    x_filtered = base_function['xaxis'].values
    y_filtered = base_function['yaxis'].values
    x_averaged = []
    y_averaged = []
    for i in np.arange(np.min(x_filtered), np.max(x_filtered), window_size):
        # get the indices of the points within the current window
        indices = np.where((x_filtered >= i) & (x_filtered < i + window_size))[0]
        if len(indices) == 0:
            continue
        # calculate the average x and y values for the points in the current window
        x_avg = np.mean(x_filtered[indices])
        y_avg = np.mean(y_filtered[indices])
        x_averaged.append(x_avg)
        y_averaged.append(y_avg)
    return x_averaged, y_averaged


# 1. Get starting base function pchip
starting_base = iteration_1['base']
# slices_base = np.arange(1700, 1850, step=1)
slices_base = np.arange(1866, 2036, step=1)
starting_base = filter_base(starting_base, slices_base)
x_averaged, y_averaged = average_base(starting_base, 0.05)
pchip_base = PchipInterpolator(x_averaged, y_averaged)

# 2. Get starting minimized_df, shifted axis and initial y values
starting_minimized = iteration_1['minimized_df'].set_index('slice').copy()
starting_axis = iteration_1['shifted'].copy()
iterations[0] = iteration_1
starting_vals = vals_it_0.copy()

# 5. Plot data from iteration 0
p1 = figure(title=f'a. Rotation angle vs. x0)', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
p2 = figure(title='c. Rotation angle vs. RMSE', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
p3 = figure(title='b. x0 vs. difference (x0-angle)', x_range=(-15, 15), y_range=(-1.5, 1.5), tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
p4 = figure(title='d. Rotation angle vs. amplitude', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
p6 = figure(title= 'Modified slices', width = 1000, height = 350, tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
p1.line(x=starting_minimized['x0'], y=starting_minimized['angle'], color=new_colors[0], legend_label=f'iteration {0}', line_width=2)
p2.line(x=starting_minimized['angle'], y=starting_minimized['rmse'], color=new_colors[0], legend_label=f'iteration {0}', line_width=2)
xaxis, ypredictions, model = polyfit(3, starting_minimized)
p3.line(x=xaxis, y=ypredictions, color=new_colors[0], legend_label=f'iteration {0}', line_width=2, dash='dashed')
p3.line(x=starting_minimized['x0'], y=starting_minimized['difference'], line_width=2, color=new_colors[0], legend_label=f'iteration {0}', )
p4.line(x=starting_minimized['angle'], y=starting_minimized['amplitude'], line_width=2, color=new_colors[0], legend_label=f'iteration {0}',)

# Plot with shifted axis
shifted_plot = figure(title='Shifted optosurf axis', width=1100, height=350, x_range=(-22,22), tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
shifted_plot.circle(x=starting_axis, y=np.zeros(len(onaxis)), line_width=2, 
                          color = new_colors[0], size = 7, legend_label=f'Iteration 0',)

# Base function plot
base_functions = figure(title = 'Base functions per iteration', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")], width = 800, height = 600, x_range = (-15,15))
base_functions.line(x=x_averaged, y=y_averaged, line_width=2, color=new_colors[0], legend_label=f'iteration 0 average',)
base_functions.triangle(x=x_averaged, y=y_averaged, size=8, color=new_colors[0], legend_label=f'iteration 0 average points',)
base_functions.circle(x=starting_base['xaxis'], y=starting_base['yaxis'], size=5, color=new_colors[0], legend_label=f'iteration {0} all points',)

shifted_axis_array = []
polynomials = {0:3, 1:3, 2:3, 3:3}
slices_amp = np.arange(slice1, slice2+1, step=base_step)
# 6. Plot the slices for iteration 0 (previously calculated)
for k, slice in enumerate(slices_amp):
    if slice % 250 == 0:
        y_exp = experimental_data[:,slice]
        y_0 = vals_it_0[:,slice]
        p6.circle(x=starting_axis, y=y_exp, size=4, color = new_colors[0], legend_label = f'Original data', )
        p6.line(x=starting_axis, y=y_exp, line_width=2, color = new_colors[0], legend_label = f'Original data')
        p6.circle(x=starting_axis, y=y_0, size=4, color = new_colors[1], legend_label = f'Iteration 0', )
        p6.line(x=starting_axis, y=y_0, line_width=2, color = new_colors[1], legend_label = f'Iteration 0')

# print('starting values max iteration 0')
# print(np.max(starting_vals))
    
# 8. Iterate minimization process with the starting base function, minimized_df and shifted axis
for u in range(1, 10):
    # 7. Calculate the p chip with the starting base function
    x_averaged, y_averaged = average_base(starting_base, 0.05)
    pchip_base = PchipInterpolator(x_averaged, y_averaged)
    
    # 8. Call the minimization procedure with the starting minimized_df and shifted axis
    min_df = minimize_procedure(slices_amp, starting_minimized, starting_axis, pchip_base, experimental_data)
    # min_df = minimize_procedure(slices_amp, starting_minimized, starting_axis, pchip_base, starting_vals)
    
    # 9. Plot rotation angle vs x0 and rotation angle vs RMSE
    p1.line(x=min_df['angle'], y=min_df['x0'], color=new_colors[u], legend_label=f'iteration {u}', line_width=2)
    p2.line(x=min_df['angle'], y=min_df['rmse'], color=new_colors[u], legend_label=f'iteration {u}', line_width=2)

    # 10. x0 vs difference polynomial fit
    poly = PolynomialFeatures(degree=1)
    x0 = min_df['angle'].values
    difference = min_df['difference'].values
    x0_poly = poly.fit_transform(x0.reshape(-1,1))
    model = LinearRegression().fit(X=x0_poly, y=difference)
    xaxis = np.arange(-15.5, 15.5, 0.1)
    xaxis_poly = poly.transform(xaxis.reshape(-1,1))
    ypredictions = model.predict(xaxis_poly)
    
    # 12. x0 vs difference plot
    p3.line(x=min_df['angle'], y=min_df['difference'], line_width=2, color=new_colors[u], legend_label=f'iteration {u}', )
    p3.line(x=xaxis, y=ypredictions, line_width=2, color = new_colors[u], dash='dashed')
    
    # 13. Estimate new shifted axis
    onaxis_poly_temp = poly.transform(starting_axis.reshape(-1,1))
    y_shifted_temp = model.predict(onaxis_poly_temp)
    p3.circle(x=starting_axis, y=y_shifted_temp, size=8, fill_color = new_colors[u], color='blue', legend_label=f'iteration {u} points',)
    shifted_axis_temp = starting_axis + y_shifted_temp
    shifted_plot.circle(x=shifted_axis_temp, y=np.zeros(len(onaxis))+0.1*u, line_width=2, 
                        color = new_colors[u], size = 7, legend_label=f'iteration {u}',)

    # 10. x0 vs difference pchip fit
    # x0 = min_df['angle'].values
    # difference = min_df['difference'].values
    # x0 = np.insert(x0, 0, -20)
    # x0 = np.append(x0, 20)
    # difference = np.insert(difference, 0, difference[0])
    # difference = np.append(difference, difference[-1])
    # pchip_diff = PchipInterpolator(x0, difference)
    # ydifferences = pchip_diff(xaxis)

    # # 12. x0 vs difference plot
    # p3.line(x=min_df['angle'], y=min_df['difference'], line_width=2, color=new_colors[u], legend_label=f'iteration {u}', )
    # # p3.line(x=xaxis, y=ydifferences, line_width=2, color = new_colors[u], dash='dashed')
    
    # # 13. Estimate new shifted axis
    # y_shifted_temp = pchip_diff(starting_axis)
    # p3.circle(x=starting_axis, y=y_shifted_temp, size=8, fill_color = new_colors[u], color='blue', legend_label=f'iteration {u} points',)
    # shifted_axis_temp = starting_axis + y_shifted_temp
    # shifted_plot.circle(x=shifted_axis_temp, y=np.zeros(len(onaxis))+0.1*u, line_width=2, 
    #                     color = new_colors[u], size = 7, legend_label=f'iteration {u}',)
   

    # 14. mechanical angle vs amplitude pchip - amplitude correction factor
    angles = min_df['angle'].values
    amplitudes = min_df['amplitude'].values
    p4.line(x=angles, y=amplitudes, line_width=2, color=new_colors[u], legend_label=f'iteration {u}',)
    angles = np.insert(angles, 0, -20)
    angles = np.append(angles, 20)
    amplitudes = np.insert(amplitudes, 0, amplitudes[0])
    amplitudes = np.append(amplitudes, amplitudes[-1])
    pchip_amp = PchipInterpolator(angles, amplitudes)
    amp =  pchip_amp(shifted_axis_temp)
    p4.circle(x=shifted_axis_temp, y=amp, size=8, color = 'blue', fill_color = new_colors[u], legend_label = f'iteration {u} pchip')

    # 15. Calculate new base function
    min_df_2 = min_df.copy().set_index('slice')
    x_base = []
    y_base = []
    slice_array = []
    temp_vals = np.zeros((32, 4096))
    # print('creating temp vals')
    # print(np.max(temp_vals))
    for k, slice in enumerate(slices_amp):
        # y_temp = starting_vals[:,slice]
        y_temp = experimental_data[:, slice]
        angle = min_df_2.loc[slice, 'angle']
        # amp = min_df_2.loc[slice, 'amplitude']
        # yamp = y_temp / amp
        yamp = y_temp*amp
        # yamp = y_temp
        if slice % 250 == 0:
            p6.circle(x=shifted_axis_temp, y=yamp, size=7, color = new_colors[u], legend_label = f'iteration {u}', )
            p6.line(x=shifted_axis_temp, y=yamp, line_width=2, color = new_colors[u], legend_label = f'iteration {u}', line_dash = 'dashed')
        
        shifted_axis_2 = shifted_axis_temp + angle
        x_base.extend(list(shifted_axis_2))
        y_base.extend(list(yamp))
        slice_array.extend(np.ones(len(shifted_axis_2))*slice)
        temp_vals[:,slice] = yamp
        # vals_new[:,slice] = yamp
    # print(f'starting values max iteration {u}')
    # print(np.max(temp_vals))
    starting_vals = temp_vals.copy()

    base_function_df = pd.DataFrame({'xaxis': x_base, 'yaxis': y_base, 'slice': slice_array})
    base_function_df = base_function_df.sort_values(by='xaxis')
    base_function_merged = base_function_df.merge(min_df, on=['slice', 'slice'])
    base_function_merged.drop(columns=['rmse', 'difference'], inplace=True)
    base_function_merged['yaxisamp'] = base_function_merged['yaxis'] / base_function_merged['amplitude']

    # 13. Update the starting axis and minimized_df
    starting_axis = shifted_axis_temp
    starting_minimized = min_df.copy().set_index('slice')
    
    sliceb = np.arange(1866, 2036, step=1)
    temp_base_function = base_function_merged.copy()
    temp_base_function = filter_base(temp_base_function, sliceb)

    base_functions.line(x=x_averaged, y=y_averaged, line_width=2, color=new_colors[u], legend_label=f'iteration {u} average',)
    base_functions.circle(x=x_averaged, y=y_averaged, size=5, color=new_colors[u], legend_label=f'iteration {u} average points',)
    base_functions.circle(x=temp_base_function['xaxis'], y=temp_base_function['yaxis'], size=5, color=new_colors[u], legend_label=f'iteration {u} all points',)
    starting_base = temp_base_function.copy()
    
    # 14. Save dictionary
    temp_dict = {'minimized_df': min_df, 'shifted': shifted_axis_temp, 'base': base_function_merged.copy()}
    iterations[u] = temp_dict

p1 = plot_format(p1, "Rotation angle", "x0", "top_left", "10pt", "11pt", "8pt")
p2 = plot_format(p2, "Rotation angle", "RMSE", "top_left", "10pt", "11pt", "8pt")
p3 = plot_format(p3, "x0", "difference", "top_left", "10pt", "11pt", "8pt")
p4 = plot_format(p4, "Rotation angle", "Amplitude", "top_left", "10pt", "11pt", "8pt")
p6 = plot_format(p6, "Angle", "Amplitude", "top_right", "10pt", "11pt", "8pt")
shifted_plot = plot_format(shifted_plot, "Angle", "", "top_left", "10pt", "11pt", "8pt")
shifted_plot.xaxis.ticker.desired_num_ticks = 40
base_functions = plot_format(base_functions, "Angle", "Amplitude", "top_left", "10pt", "11pt", "8pt")
grid = gridplot(children=[p1, p2, p3, p4], ncols=2, merge_tools=False, width = 650, height = 500)
show(grid)
show(p6)
Code
show(shifted_plot)
Code
show(base_functions)
Code
axis = iterations[5]['shifted']
base = iterations[5]['base']
base_functions_starting = {}
base_plots = []
colorn = 0

step = 170
# for j in range(1182, 2700, step):
for j in range(1100, 2700, step):
    base_plot = figure(title = f'{offaxis[j]:.3f} deg ({j}) to {offaxis[j+step]:.3f} deg ({j+step})', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")], 
                       width = 650, height = 500, x_range = (-5.0, 5.0), y_range = (0, 40000))
    
    # 1 degree sections
    filter_section = np.arange(j, j+step, step=1)
    base_sorted = base[base['slice'].isin(filter_section)].sort_values(by=['xaxis']).reset_index().copy()
    base_sorted.loc[0] = [0, -500, 0, 0, 0, 0, 0, 0]
    base_sorted.loc[-1] = [0, 500, 0, 0, 0, 0, 0, 0]
    
    # # average window
    x_filtered = base_sorted['xaxis'].values 
    y_filtered = base_sorted['yaxis'].values
    window_size = 0.05
    x_averaged = []
    y_averaged = []
    for i in np.arange(np.min(x_filtered), np.max(x_filtered), window_size):
        # get the indices of the points within the current window
        indices = np.where((x_filtered >= i) & (x_filtered < i + window_size))[0]
        if len(indices) == 0:
            continue
        # calculate the average x and y values for the points in the current window
        x_avg = np.mean(x_filtered[indices])
        y_avg = np.mean(y_filtered[indices])
        x_averaged.append(x_avg)
        y_averaged.append(y_avg)
    x_averaged = np.array(x_averaged)
    y_averaged = np.array(y_averaged)
    pchip = PchipInterpolator(x_averaged, y_averaged)
    # cubic = CubicSpline(x_averaged, y_averaged)
    x_interpolated = np.arange(-20, 20, step=0.01)
    y_interpolated = pchip(x_interpolated)
    # y_interpolated = cubic(x_interpolated)
    base_plot.line(x= x_interpolated, y=y_interpolated, line_width=7, legend_label = 'cubic line',)
    # base_plot.line(x= np.arange(-20, 20, step=0.1), y=pchip(np.arange(-20, 20, step=0.1)), line_width=7, legend_label = 'Pchip line',)
    base_plot.circle(x= x_averaged, y=y_averaged, size=7, legend_label = 'Average points', color = 'red')
    base_functions_starting[colorn] = [x_averaged, y_averaged]
    
    # pchip = PchipInterpolator(base_sorted['xaxis'], base_sorted['yaxisamp'])

    # Plot for two 0.5 degree sections
    filter1 = np.arange(j, j+step//2, step=1)
    filter2 = np.arange(j+step//2, j+step, step=1)
    filters = [filter1, filter2]
    for q, filter in enumerate(filters):

        base_sorted = base[base['slice'].isin(filter)]
        base_sorted = base_sorted.sort_values(by=['xaxis'])
        # base_plot.line(x=base_sorted['xaxis'], y=base_sorted['yaxisamp'], line_width=5, 
        #              color=new_colors[colorn+q],)
        base_plot.circle(x=base_sorted['xaxis'], y=base_sorted['yaxis'], size=7, 
                     color=new_colors[colorn+q], legend_label = f'{offaxis[filter[0]]:.3f} to {offaxis[filter[-1]]:.3f} deg',)
        base_plot = plot_format(base_plot, "Angle", "Amplitude", "top_right", "10pt", "11pt", "8pt")
    base_plots.append(base_plot)
    colorn += 1

grid = gridplot(children=base_plots, ncols=3, merge_tools=False, width = 580, height = 400)
show(grid)
Code
all_filtered = np.arange(1100, 2700, step=1)
base_section = base[base['slice'].isin(all_filtered)]
base_section = base_section.sort_values(by=['xaxis'])
# max = base_section['yaxisamp'].max()

filter = np.arange(1100, 2700, step=1)
base_plot_filtered = figure(title = f'{offaxis[filter[0]]:.3f} to {offaxis[filter[-1]]:.3f} deg', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")], 
                       width = 580, height = 400, x_range = (-20.0, 20.0),)
base_sorted = base[base['slice'].isin(filter)]
base_sorted = base_sorted.sort_values(by=['xaxis'])
base_plot_filtered.circle(x=base_sorted['xaxis'], y=base_sorted['yaxisamp'], size=5, 
                     color=new_colors[0],)
base_plot_filtered = plot_format(base_plot_filtered, "Angle", "Amplitude", "top_right", "10pt", "11pt", "8pt")
show(base_plot_filtered)
Code
axis = iterations[5]['shifted']
base = iterations[5]['base']

def moving_average(x, window_size):
    weights = np.repeat(1.0, window_size) / window_size
    return np.convolve(x, weights, mode='valid')

x_smooth = moving_average(base_sorted['xaxis'], 5)
y_smooth = moving_average(base_sorted['yaxis'], 5)

filter_section = np.arange(1950, 2120, step=1)
base_sorted = base[base['slice'].isin(filter_section)].sort_values(by=['xaxis']).reset_index().copy()
base_sorted.loc[0] = [0, -500, 0, 0, 0, 0, 0, 0]
base_sorted.loc[-1] = [0, 500, 0, 0, 0, 0, 0, 0]
# x_base = base_sorted['xaxis']
# y_base = base_sorted['yaxis']
x_base, y_base = average_base(base_sorted, 0.05)
pchipbase = PchipInterpolator(x_base, y_base)
x_interpolated = np.arange(-16.5, 16.5, step=0.01)
y_interpolated = pchipbase(x_interpolated)


base_plot = figure(title = f'{offaxis[filter_section[0]]:.3f} to {offaxis[filter_section[-1]]:.3f} deg', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")], 
                   x_range = (-3.5, 3.5), y_range = (0, 40000), width = 580, height = 400)
base_plot.line(x= x_interpolated, y=y_interpolated, line_width=7, legend_label = 'pchip', color = new_colors[5])
base_plot.circle(x=base_sorted['xaxis'], y=base_sorted['yaxis'], size=5, color = new_colors[2], legend_label = 'Base function')
base_plot.circle(x=x_base, y=y_base, size=5, color = new_colors[3], legend_label = 'Average points')


base_plot = plot_format(base_plot, "Angle", "Amplitude", "top_left", "10pt", "11pt", "8pt")
show(base_plot)
Code
path2 = '/Users/carlosreyes/Desktop/Academics/postdoc/confirm/python/optosurf/Streamlit/data/f/'
file1 = 'Rotate_Ann7_offaxis_10degscan.csv'
file2 = 'Rotate_Ann7_onaxis_10degscan.csv' 
files = [file1, file2,]

offaxis_ann7 = np.arange(-512, 512, 1)*20/1024 + 1.82
onaxis_ann7  = iterations[4]['shifted']

subplot, spots_ann7 = subplot3d(files, offaxis, onaxis, 1, 2, path2)
subplot.update_scenes(camera_eye=dict(x=-0.75, y=-1.0, z=1.25))
subplot.update_layout(width=1200, height=550)
subplot
Unable to display output for mime type(s): application/vnd.plotly.v1+json
Code
def cost_function(params, y, x, phcip):
    """
    Generates a 3D plot of the data

    Parameters
    ----------
    params (list): List of parameters to fit
    y (array): Array with the data to fit
    
    Returns
    -------
    rmse (float): Root mean square error of the fit0.
    """
    x0, A0 = params
    x_new = x + x0
    y_modified = A0*pchip(x_new)

    mse = np.mean((y - y_modified) ** 2)
    rmse = np.sqrt(mse)
    return rmse

slices_ann7 = np.arange(10, 1000, step=1)
minimized_df_ann7 = pd.DataFrame(columns=["angle", "x0", "difference", "slice", "amplitude", "rmse"])
pchipbase = PchipInterpolator(x_base, y_base)

offaxis_ann7 = np.arange(-512, 512, 1)*20/1024 + 1.82
onaxis_ann7  = iterations[4]['shifted']
data_ann7 = spots_ann7['Rotate_Ann7_onaxis_10degscan.csv']

for j, slice in enumerate(slices_ann7): 
        x0 = offaxis_ann7[slice]
        Abase = 1.0
        guess_no_back = [x0, Abase]

        # Call minimization function
        y = data_ann7[:,slice]
        cost_fn = lambda p:cost_function(p, y, onaxis_ann7, pchipbase)

        result = minimize(cost_fn, guess_no_back, method=method,)
        optimized_parameters = list(result.x)
        x0_opt, A0_opt, = optimized_parameters

        # Calculate optimized function
        x_new_opt = onaxis_ann7 + (x0_opt)
        # y_optimized = A0_opt*pchip(x_new_opt) 
        y_optimized = A0_opt*pchip(x_new_opt)

        # Calculate error
        mse = np.mean((y - y_optimized) ** 2)
        rmse = np.sqrt(mse)
        row = [offaxis_ann7[slice], x0_opt, offaxis_ann7[slice] - x0_opt, slice, A0_opt, rmse]
        minimized_df_ann7.loc[j] = row
Code
min_df_iteration5 = iterations[4]['minimized_df']
p1_ann7 = figure(title=f'a. Rotation angle vs. x0 (ann7)')
p2_ann7 = figure(title='c. Rotation angle vs. RMSE (ann7)')
p3_ann7 = figure(title='b. x0 vs. difference (mechanical angle - x0) (ann7)')
p4_ann7 = figure(title='d. Rotation angle vs. amplitude (ann7)')
p5_ann7 = figure(title='Shifted optosurf axis (ann7)', width=1000, height=350, x_range=(-25,25))

def linear_function(x, m , b):
    return m*x + b
shifted_axis_dict = {}
pchip_amps = {}
angles_interp_dict = {}
amplitudes_interp_dict = {}


# a. Rotation angle vs. x0
p1_ann7.line(x=min_df_iteration5['angle'], y=min_df_iteration5['x0'], line_width=2, color = new_colors[0], legend_label = 'ann10 reference')
p1_ann7.line(x=minimized_df_ann7['angle'], y=minimized_df_ann7['x0'], line_width=2, color = new_colors[1], legend_label = 'ann7 - iteration 0')

# b. Rotation angle vs RMSE
p2_ann7.line(x=min_df_iteration5['angle'], y=min_df_iteration5['rmse'], line_width=2, color = new_colors[0], legend_label = 'ann10 reference' )
p2_ann7.line(x=minimized_df_ann7['angle'], y=minimized_df_ann7['rmse'], line_width=2, color = new_colors[1], legend_label = 'ann7 - iteration 0')

# c. Polynomial fit of x0 vs difference
x0 = min_df_iteration5['x0'].values
angles = min_df_iteration5['angle'].values
difference = min_df_iteration5['difference'].values
p3_ann7.line(x=angles, y=difference, line_width=2, color = new_colors[0], legend_label = 'ann10 reference' )

x0_ann7 = minimized_df_ann7['x0'].values
angles_ann7 = minimized_df_ann7['angle'].values
difference_ann7 = minimized_df_ann7['difference'].values
p3_ann7.line(x=angles_ann7, y=difference_ann7, line_width=2, color = new_colors[1], legend_label = 'ann7 - iteration 0')

# e. Rotation angle vs amplitude with interpolation
angles = min_df_iteration5['angle'].values
amplitudes = min_df_iteration5['amplitude'].values
p4_ann7.line(x=angles, y=amplitudes, line_width=1, color = new_colors[0], legend_label = 'ann10 reference')

angles_ann7 = minimized_df_ann7['angle'].values
amplitudes_ann7 = minimized_df_ann7['amplitude'].values
p4_ann7.line(x=angles_ann7, y=amplitudes_ann7, line_width=1, color = new_colors[1], legend_label = 'ann7 - iteration 0')


# f. Plots format
p1_ann7.xaxis.ticker.desired_num_ticks = 10
p1_ann7 = plot_format(p1_ann7, "Rotation angle", "x0", "top_left", "10pt", "8pt", "8pt")
p2_ann7.xaxis.ticker.desired_num_ticks = 10
p2_ann7 = plot_format(p2_ann7, "Rotation angle", "RMSE", "top_left", "10pt", "11pt", "8pt")
p3_ann7.xaxis.ticker.desired_num_ticks = 10
p3_ann7 = plot_format(p3_ann7, "Rotation angle", "difference", "top_left", "10pt", "11pt", "8pt")
p4_ann7.xaxis.ticker.desired_num_ticks = 10
p4_ann7 = plot_format(p4_ann7, "Rotation angle", "amplitude", "top_left", "10pt", "11pt", "8pt")


grid_ann7 = gridplot(children=[p1_ann7, p2_ann7, p3_ann7, p4_ann7], ncols=2, merge_tools=False, width = 550, height = 400)
show(grid_ann7)

4 Rough data optimization

Code
from bokeh.palettes import Set3
from IPython.display import HTML

# Get base function
x_base_a = x_base
y_base_a = y_base
pchip = pchipbase
shifted_axis = iterations[5]['shifted']

# shifted_axis = onaxis
# Get rough data
rough_df = pd.read_excel('data/rough_samples_shifted.xlsx', sheet_name='Data')
# rough_df['shiftedaxis'] = shifted_axis
# source_rough = ColumnDataSource(rough_df)
x_rough = rough_df['xaxis'].values.round(3)

# Get initial guess
guess_df = pd.read_excel('data/rough_samples_shifted.xlsx', sheet_name='SuperGaussian')
guess_df = guess_df.set_index('Variables')

# Create df that will save optmized parameters
columns = ['ann1_opt', 'pt2_opt', 'pt2b_opt', 'pt2c_opt', 'pt2d_opt', 'pt2e_opt']
index = ['x0', 'Abase', 'sigma', 'Agaussian', 'n']
optimized_df = pd.DataFrame(columns=columns, index=index)

# Define gaussian function
supergaussian = lambda x, x0, sigma, A1, n: A1 * np.exp(-abs(((x-x0)/sigma))**n)

# 6. Define cost function
def cost_function(params, y):
    x0, A0, sigma, A1, n = params
    # Get new x axis
    x_new = shifted_axis + x0 
    # interpolate base function with respect to x_new (32 points)
    y_base_modified = A0*pchip(x_new) 
    # calculate background on original axis and with x0
    y_background = supergaussian(x_new, x0, sigma, A1, n)
    # calculate modified function
    y_modified = y_base_modified + y_background

    mse = np.mean((y - y_modified) ** 2)
    rmse = np.sqrt(mse)
    # Compare directly with 32 points experimental data
    return rmse

# 7. Iterate over the experimental data
rough_plots = []
color_palette = Set3[len(rough_df.columns[1:])+2]
bounds = ((None, None), (0, None), (0, None), (0, None))
backgrounds = figure(title = 'Background functions with experimental data', width = 550, height = 450, tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
backgroundsbg = figure(title = 'Background functions', width = 550, height = 450, tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
downsamplesg = figure(title = 'Gaussian downsampled points', width = 550, height = 450, tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
differenceg = figure(title = 'Gaussian Experimental vs Optimized differences', width = 550, height = 450, tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])

for i, col in enumerate(rough_df.columns[1:]):
    # 8. Get initial guesses
    rough_plot = figure(title = str(col), width = 550, height = 450, tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
    guess = [guess_df.loc[var][col] for var in ['x0', 'Abase', 'sigma', 'Agaussian', 'n']]

    # 8. Call minimization function
    y_rough = rough_df[col].copy().values
    # y_rough = y_rough
    cost_fn = lambda p:cost_function(p, y_rough)
    result = minimize(cost_fn, guess)
    # result = minimize(cost_fn, params, bounds=bounds)
    optimized_parameters = result.x
    colu = col + '_opt'
    x0_opt, A0_opt, sigma_opt, A1_opt, n_opt = optimized_parameters
    optimized_df.loc['x0'][colu] = x0_opt
    optimized_df.loc['Abase'][colu] = A0_opt
    optimized_df.loc['sigma'][colu] = sigma_opt
    optimized_df.loc['Agaussian'][colu] = A1_opt
    optimized_df.loc['n'][colu] = n_opt

    # 7. Calculate new optimized modified function
    x_new_opt = shifted_axis + x0_opt 
    # interpolate base function with respect to x_new (32 points)
    y_base_opt = A0_opt*pchip(x_new_opt) 
    # calculate background on original axis and with x0
    y_background_opt = supergaussian(x_new_opt, x0_opt, sigma_opt, A1_opt, n_opt)
    # calculate optmized function
    y_optimized = y_base_opt + y_background_opt

    vline = Span(location=0.0, dimension = 'height', line_color='#FEEED9', line_width=1)
    rough_plot.add_layout(vline)

    # Plot optimize function lines
    rough_plot.line(shifted_axis, y_base_opt, legend_label = 'Base', line_width = 5, color='#F96F5D')
    rough_plot.line(shifted_axis, y_background_opt, legend_label = 'Bbackground', line_width = 5, color='#F9B5AC')
    rough_plot.line(shifted_axis, y_optimized, legend_label = 'Optimized function', line_width = 5, color='#987284')
    rough_plot.triangle(shifted_axis, y_optimized, legend_label = 'Optimized points', size = 8, color=color_palette[1])
    backgrounds.line(shifted_axis, y_background_opt, color = color_palette[i], line_width = 5 , legend_label = f"Background {col}")
    backgroundsbg.line(shifted_axis, y_background_opt, color = color_palette[i], line_width = 5 , legend_label = f"Background {col}")
    downsamplesg.line(shifted_axis, y_optimized, line_width=4, legend_label = f'Downsampling {col}', color = color_palette[i+1],  alpha = 0.9, line_dash='dashed')
    downsamplesg.triangle(shifted_axis, y_optimized, size = 13, legend_label = f'Downsampling {col}', color = color_palette[i+1])

    # Plot rough experimental data
    rough_plot.line(shifted_axis, rough_df[col], color = '#9DC3E6', legend_label = str(col), line_width=4, line_dash = 'dashed')
    rough_plot.circle(shifted_axis, rough_df[col], fill_color= color_palette[i], size=7, legend_label = f"{col} points")
    backgrounds.line(shifted_axis, rough_df[col], color = color_palette[i], legend_label = str(col), line_width=4)
    backgrounds.circle(shifted_axis, rough_df[col], fill_color= color_palette[i], size=7, legend_label = str(col))
    downsamplesg.line(shifted_axis, rough_df[col], color = color_palette[i], legend_label = str(col), line_width=4)
    downsamplesg.circle(shifted_axis, rough_df[col], fill_color= color_palette[i], size=7, legend_label = str(col))
    
    # Plot format
    # rough_plot.y_range = Range1d(-5000, 50000)
    rough_plot.xaxis.ticker.desired_num_ticks = 10
    rough_plot.yaxis.ticker.desired_num_ticks = 10
    rough_plot = plot_format(rough_plot, "Degrees", "Intensity", "top_left", "10pt", "8pt", "9pt")
    rough_plots.append(rough_plot)

    # Difference plot
    diff = y_rough - y_optimized
    differenceg.line(x=x_rough, y=diff, legend_label = col, color = color_palette[i], line_width=4)
    differenceg.circle(x=x_rough, y=diff, legend_label = col, fill_color= color_palette[i], size=7)

plots = [backgrounds, backgroundsbg, downsamplesg, differenceg]
for plot in plots:
    plot = plot_format(plot, "Degrees", "Intensity", "top_left", "10pt", "10pt", "9pt")
    # plot.y_range = Range1d(-5000, 50000)
    rough_plots.append(plot)
    plot.xaxis.ticker.desired_num_ticks = 10
    plot.yaxis.ticker.desired_num_ticks = 10
# backgroundsbg.y_range = Range1d(-1000, 6000)
# differenceg.y_range = Range1d(-3000, 3000)
downsamplesg.add_layout(vline)
backgrounds.add_layout(vline)
backgroundsbg.add_layout(vline)

merged_df_2 = guess_df.join(optimized_df)\
    [['ann1', 'ann1_opt', 'pt2', 'pt2_opt', 'pt2b', 'pt2b_opt', 'pt2c', 'pt2c_opt', 'pt2d', 'pt2d_opt', 'pt2e', 'pt2e_opt']]

merged_df_html_2 = merged_df_2.style.set_table_attributes('style="font-size: 11px"').render()
display(HTML(merged_df_html_2))
# base_function = figure(title = f'base function', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")], 
#                        width = 650, height = 500, x_range = (-20.0, 20.0),)
# base_function.line(x= np.arange(-20, 20, step=0.1), y=pchip(np.arange(-20, 20, step=0.1)), line_width=7, legend_label = 'Pchip line',)
# base_function.circle(x=x_base_a, y=y_base_a, size=5, color=new_colors[2])

# base_function = plot_format(base_function, "Angle", "Amplitude", "top_right", "10pt", "11pt", "8pt")
# show(base_function)
KeyError: 5
Code
grid_rough = gridplot(children = rough_plots, ncols = 3, merge_tools=False, width = 480, height = 430)
show(grid_rough)
Code
# ranges = [[1690, 1774], [1775, 1859], [1860, 1944], [1945, 2029]]
Code
print(base_sorted.info)
<bound method DataFrame.info of        index       xaxis     yaxis   slice     angle        x0  amplitude  \
 0         0 -500.000000  0.000000     0.0  0.000000  0.000000   0.000000   
 1     59232  -15.848137  1.015369  1951.0  0.491641  0.491433   0.999717   
 2     59264  -15.842278  2.030738  1952.0  0.497500  0.497043   0.999773   
 3     59296  -15.836418  4.061475  1953.0  0.503359  0.502753   0.999934   
 4     59328  -15.830559  1.015369  1954.0  0.509219  0.509535   1.000406   
...      ...         ...       ...     ...       ...       ...        ...   
 5436  64543   17.767236  0.995634  2116.0  1.458438  1.461348   0.999079   
 5437  64575   17.773095  0.000000  2117.0  1.464297  1.466810   0.999198   
 5438  64607   17.778955  1.991267  2118.0  1.470156  1.471987   0.999179   
 5439  64639   17.784814  2.986901  2119.0  1.476016  1.476681   0.999793   
-1         0  500.000000  0.000000     0.0  0.000000  0.000000   0.000000   

       yaxisamp  
 0     0.000000  
 1     1.015656  
 2     2.031199  
 3     4.061744  
 4     1.014956  
...         ...  
 5436  0.996552  
 5437  0.000000  
 5438  1.992904  
 5439  2.987518  
-1     0.000000  

[5441 rows x 8 columns]>
Code
print((offaxis[1] - offaxis[0]))
print(0.5//(offaxis[1] - offaxis[0]))
0.005859375
85.0
Code
print((offaxis[2050]))
1.07171875
Code
# #| column: screen
# iteration_plots = []
# for iteration in iterations:
#     base_functions_plot_2 = figure(title=f'Base functions iteration {iteration}', tooltips = [("index", "$index"),("(x,y)", "($x, $y)")])
#     base_temp = iterations[iteration]['base']
#     base_functions_plot_2.circle(x=base_temp['xaxis'], y=base_temp['yaxisamp'], size=5, color=new_colors[iteration],)
#     iteration_plots.append(base_functions_plot_2)
#     base_functions_plot_2 = plot_format(base_functions_plot_2, "Angle", "Amplitude", "top_left", "10pt", "11pt", "8pt")
# grid_base = gridplot(children=iteration_plots, ncols=3, merge_tools=False, width = 550, height = 400)
# show(grid_base)
Code
vals = spots['refrockley_focus_0.csv']

for k, slice in enumerate(slices_base):
    # some code
    vals[:,slice] = yamp